home *** CD-ROM | disk | FTP | other *** search
/ Amiga Format CD 24 / Amiga Format AFCD24 (Feb 1998, Issue 108).iso / -seriously_amiga- / shareware / programming / c / amivogl-mdev / src / patches.c < prev    next >
C/C++ Source or Header  |  1998-01-12  |  11KB  |  541 lines
  1. #include "vogl.h"
  2.  
  3. static int    ntsegs;
  4. static int  tsegs = 10;
  5. static int  nusegs;
  6. static int  usegs = 10;
  7. static int  ntcurves = 10;
  8. static int  nucurves = 10;
  9. static int  ntiter = 1; 
  10. static int  nuiter = 1;
  11. static int  makeprec_called = 0;
  12.  
  13. static Matrix et;
  14. static Matrix eu;
  15. static Matrix ones =    {
  16.     {1.0, 1.0, 1.0, 1.0}
  17.     ,
  18.     {1.0, 1.0, 1.0, 1.0}
  19.     ,
  20.     {1.0, 1.0, 1.0, 1.0}
  21.     ,
  22.     {1.0, 1.0, 1.0, 1.0}
  23.     };
  24.  
  25. /* ---------------------------------------------------------------------
  26.  * Prototypes:
  27.  */
  28. #ifdef __PROTOTYPE__
  29.  
  30. static void makeprec(void);                            /* patches.c       */
  31. static void replace( Tensor, Matrix, int);             /* patches.c       */
  32. static void iterate( register Tensor, int);            /* patches.c       */
  33. static void extract( register Matrix, register Tensor, /* patches.c       */
  34.    register int);
  35. static float addemup(Matrix);                          /* patches.c       */
  36.  
  37. #else    /* __PROTOTYPE__ */
  38.  
  39. static void makeprec();                                /* patches.c       */
  40. static void replace();                                 /* patches.c       */
  41. static void iterate();                                 /* patches.c       */
  42. static void extract();                                 /* patches.c       */
  43. static float addemup();                                /* patches.c       */
  44.  
  45. #endif    /* __PROTOTYPE__ */
  46.  
  47.  
  48. /*
  49.  * defbasis
  50.  *
  51.  *    Specify a basis matrix for a patch or curve.
  52.  */
  53. void defbasis(
  54.   short id,
  55.   Matrix mat)
  56. {
  57.         if(!vdevice.initialised)
  58.                 verror("defbasis: vogl not initialised");
  59.  
  60.         copymatrix(vdevice.bases[id], mat);
  61. }
  62.  
  63. /* ------------------------------------------------------------------------ */
  64.  
  65. /*
  66.  * patchbasis
  67.  *
  68.  *    Specify the two basis matrices for a patch
  69.  */
  70. void patchbasis(
  71.   long tb,
  72.   long ub)
  73. {
  74.         if(!vdevice.initialised)
  75.                 verror("patchbasis: vogl not initialised");
  76.  
  77.         copymatrix(vdevice.tbasis, vdevice.bases[tb]);
  78. copytranspose(vdevice.ubasis, vdevice.bases[ub]);
  79. }
  80.  
  81. /* ------------------------------------------------------------------------ */
  82.  
  83. /*
  84.  * patchcurves
  85.  *
  86.  *    Specify the number of curves to be drawn in each direction 
  87.  *    on a patch.
  88.  */
  89. void patchcurves(
  90.   long nt,
  91.   long nu)
  92. {
  93.  
  94.         if (!vdevice.initialised)
  95.                 verror("patchcurves: vogl not initialised");
  96.  
  97.         if(nt > 0 && nu > 0) {
  98.                     ntcurves = nt;
  99.                     nucurves = nu;
  100.             }
  101. else
  102.                 verror("patchcurves: number of patch curves <= 0");
  103.  
  104. /*
  105.  * Set up the difference matrices
  106.  */
  107. makeprec();
  108. }
  109.  
  110. /* ------------------------------------------------------------------------ */
  111.  
  112. /*
  113.  * patchprecision
  114.  *
  115.  *    Specify the lower limit on the number of line segments used
  116.  * to draw a curve as part of a patch. The actual number used varies
  117.  * with the number of curve segments in the "other" direction.
  118.  */
  119. void patchprecision(
  120.   long tseg,
  121.   long useg)
  122. {
  123.         if (!vdevice.initialised)
  124.                 verror("patchprecision: vogl not initialised");
  125.  
  126.         if(tseg > 0 && useg > 0) {
  127.                     tsegs = tseg;
  128.                     usegs = useg;
  129.             }
  130. else
  131.                 verror("patchprecision: number of segments <= 0");
  132. /*
  133.  * Set up the difference matrices
  134.  */
  135. makeprec();
  136. }
  137.  
  138. /* ------------------------------------------------------------------------ */
  139.  
  140. /*
  141.  * makeprec
  142.  *
  143.  *    Makes up the two precision matrices for a patch
  144.  */
  145. static void makeprec(void)
  146. {
  147. float    n2, n3;
  148.  
  149. /*
  150.  * Find ntsegs, nusegs, ntiter, nuiter....
  151.  * ie. the actual number of curve segments of a tcurve,
  152.  *     the actual number of curve segments of a ucurve,
  153.  *     and the number of times to iterate in each direction.
  154.  */
  155.  
  156. ntsegs = tsegs;
  157. ntiter = ntsegs / (nucurves - 1);
  158. if (ntsegs > ntiter * (nucurves - 1)) 
  159. ntsegs = (++ntiter) * (nucurves - 1);
  160.  
  161. nusegs = usegs;
  162. nuiter = nusegs / (ntcurves - 1);
  163. if (nusegs > nuiter * (ntcurves - 1)) 
  164. nusegs = (++nuiter) * (ntcurves - 1);
  165.  
  166. /*
  167.  * Doing the t precision matrix.....
  168.  */
  169. identmatrix(et);
  170. n2 = (float)(ntsegs * ntsegs);
  171. n3 = (float)(ntsegs * n2);
  172.  
  173. et[0][0] = et[2][2] = et[3][3] = 0.0;
  174. et[1][0] = 1.0 / n3;
  175. et[1][1] = 1.0 / n2;
  176. et[2][0] = et[3][0] = 6.0 / n3;
  177. et[2][1] = 2.0 / n2;
  178. et[1][2] = 1.0 / (float)ntsegs;
  179. et[0][3] = 1.0;
  180.  
  181. /*
  182.  * Make the Transpose of eu
  183.  */
  184. identmatrix(eu);
  185. n2 = (float)(nusegs * nusegs);
  186. n3 = (float)(nusegs * n2);
  187.  
  188. eu[0][0] = eu[2][2] = eu[3][3] = 0.0;
  189. eu[0][1] = 1.0 / n3;
  190. eu[1][1] = 1.0 / n2;
  191. eu[0][2] = eu[0][3] = 6.0 / n3;
  192. eu[1][2] = 2.0 / n2;
  193. eu[2][1] = 1.0 / (float)nusegs;
  194. eu[3][0] = 1.0;
  195.  
  196. makeprec_called = 1;
  197. }
  198.  
  199. /* ------------------------------------------------------------------------ */
  200.  
  201. /*
  202.  * patch
  203.  *
  204.  *    Draws a bicubic patch. (ie. all the w coords a 1 and the
  205.  *    basis matrices don't change that)
  206.  */
  207. void patch(
  208.   Matrix geomx,
  209.   Matrix geomy,
  210.   Matrix geomz)
  211. {
  212. rpatch(geomx, geomy, geomz, ones);
  213. }
  214.  
  215. /* ------------------------------------------------------------------------ */
  216.  
  217. /*
  218.  * rpatch
  219.  *
  220.  *    Draws rational bicubic patches.
  221.  *
  222.  *    Reference: J. H. Clark, Parametric Curves, Surfaces and volumes in 
  223.  *    computer graphics and computer aided Geometric Design.
  224.  *    Technical report No. 221, Nov 1981.
  225.  *    Computer Systems Lab. Dept's of Elecrical Eng. and Computer Science,
  226.  *    Standford University, Standford, California 94305.
  227.  */
  228. void rpatch(
  229.   Matrix geomx,
  230.   Matrix geomy,
  231.   Matrix geomz,
  232.   Matrix geomw)
  233. {
  234.  
  235. Tensor    S, R;
  236.         Matrix    tmp, tmp2;
  237. float    xlast, ylast, zlast;
  238.         int    i, j;
  239. Token    *tok;
  240.  
  241.         if (!vdevice.initialised)
  242.                 verror("patch: vogl not initialised");
  243.  
  244. /* 
  245.  *  Form S = et . tbasis . Gtensor . ubasisT . euT
  246.  */
  247.  
  248. if (!makeprec_called)
  249. makeprec();
  250.  
  251. mult4x4(tmp, et, vdevice.tbasis);
  252. mult4x4(tmp2, vdevice.ubasis, eu);
  253.  
  254. /*
  255.  * Load the geometry matrices into S.
  256.  */
  257. for (i = 0; i < 4; i++)
  258. for (j = 0; j < 4; j++) {
  259.     S[0][i][j] = geomx[i][j];
  260.     S[1][i][j] = geomy[i][j];
  261.     S[2][i][j] = geomz[i][j];
  262.     S[3][i][j] = geomw[i][j];
  263.     }
  264.  
  265. premulttensor(R, vdevice.tbasis, S);
  266. multtensor(S, vdevice.ubasis, R);
  267.  
  268. /*
  269.  * Find the last point on the curve.
  270.  */
  271. xlast = addemup(S[0]);
  272. ylast = addemup(S[1]);
  273. zlast = addemup(S[2]);
  274.  
  275. /*
  276.  * Multiply the precision matrices in.
  277.  */
  278. premulttensor(R, et, S);
  279. multtensor(S, eu, R);
  280.  
  281. if (vdevice.inobject) {
  282.     tok = newtokens(74);
  283.     tok[0].i = RPATCH;
  284.     tok[1].f = xlast;
  285.     tok[2].f = ylast;
  286.     tok[3].f = zlast;
  287.     tok[4].i = ntcurves;
  288.     tok[5].i = nucurves;
  289.     tok[6].i = ntsegs;
  290.     tok[7].i = nusegs;
  291.     tok[8].i = ntiter;
  292.     tok[9].i = nuiter;
  293.  
  294.     tok += 10;
  295.     for (i = 0; i < 4; i++)
  296.     for (j = 0; j < 4; j++) {
  297.         (tok++)->f = S[0][i][j];
  298.         (tok++)->f = S[1][i][j];
  299.         (tok++)->f = S[2][i][j];
  300.         (tok++)->f = S[3][i][j];
  301.         }
  302.  
  303.     return;
  304.     }
  305.  
  306. /*
  307.  * Multiply by the current transformation....
  308.  */
  309. transformtensor(S, vdevice.transmat->m);
  310.  
  311. /*
  312.  * Draw the patch....
  313.  */
  314. drpatch(S, ntcurves, nucurves, ntsegs, nusegs, ntiter, nuiter);
  315.  
  316. /*
  317.  * Set the current (untransformed) world spot....
  318.  */
  319. vdevice.cpW[V_X] = xlast;
  320. vdevice.cpW[V_Y] = ylast;
  321. vdevice.cpW[V_Z] = zlast;
  322. }
  323.  
  324. /* ------------------------------------------------------------------------ */
  325.  
  326. /*
  327.  * transformtensor
  328.  *
  329.  *    Transform the tensor S by the matrix m
  330.  */
  331. void transformtensor(
  332.   Tensor S,
  333.   Matrix m)
  334. {
  335. Matrix    tmp, tmp2;
  336. register    int    i;
  337.  
  338. for (i = 0; i < 4; i++) {
  339.     extract(tmp, S, i);
  340.     mult4x4(tmp2, tmp, m);
  341.     replace(S, tmp2, i);
  342.     }
  343. }
  344.  
  345. /* ------------------------------------------------------------------------ */
  346.  
  347. /*
  348.  * replace
  349.  *
  350.  *    Does the reverse of extract.
  351.  */
  352. static void replace(
  353.   Tensor a,
  354.   Matrix b,
  355.   int k)
  356. {
  357. int    i, j;
  358.  
  359. /*
  360.  * Not unwound because it only gets called once per patch.
  361.  */
  362. for (i = 0; i < 4; i++)
  363. for (j = 0; j < 4; j++)
  364. a[j][i][k] = b[i][j];
  365. }
  366.  
  367. /* ------------------------------------------------------------------------ */
  368.  
  369. /*
  370.  * drpatch
  371.  *
  372.  *    Actually does the work of drawing a patch.
  373.  */
  374. void drpatch(
  375.   Tensor R,
  376.   int ntcurves,
  377.   int nucurves,
  378.   int ntsegs,
  379.   int nusegs,
  380.   int ntiter,
  381.   int nuiter)
  382. {
  383. Tensor    S;
  384. Matrix    tmp;
  385. int    i;
  386.  
  387. /*
  388.          *  Copy R transposed into S
  389.          */
  390. copytensortrans(S, R);
  391.  
  392. for (i = 0; i < ntcurves; i++) {
  393.     extract(tmp, R, 0);
  394.     drcurve(ntsegs, tmp);
  395.     iterate(R, nuiter);
  396.     }
  397.  
  398. /*
  399.  * Now using S...
  400.  */
  401. for (i = 0; i < nucurves; i++) {
  402.     extract(tmp, S, 0);
  403.     drcurve(nusegs, tmp);
  404.     iterate(S, ntiter);
  405.     }
  406. }
  407.  
  408. /* ------------------------------------------------------------------------ */
  409.  
  410. /*
  411.  * iterate
  412.  *
  413.  *    Iterates the forward difference tensor R
  414.  */
  415. static void iterate(
  416.   register Tensor R,
  417.   int n)
  418. {
  419. register int    it;
  420.  
  421. /*
  422.  * Anyone for an unwound loop or two???
  423.  *
  424.  *    for (it = 0; it < n; it++) {
  425.  *        for (i = 0; i < 4; i++) 
  426.  *            for (j = 0; j < 4; j++)
  427.  *                for (k = 0; k < 3; k++)
  428.  *                    R[i][j][k] += R[i][j][k+1];
  429.  *    }
  430.  */
  431. for (it = 0; it < n; it++) {
  432.     R[0][0][0] += R[0][0][1];
  433.     R[0][0][1] += R[0][0][2];
  434.     R[0][0][2] += R[0][0][3];
  435.  
  436.     R[0][1][0] += R[0][1][1];
  437.     R[0][1][1] += R[0][1][2];
  438.     R[0][1][2] += R[0][1][3];
  439.  
  440.     R[0][2][0] += R[0][2][1];
  441.     R[0][2][1] += R[0][2][2];
  442.     R[0][2][2] += R[0][2][3];
  443.  
  444.     R[0][3][0] += R[0][3][1];
  445.     R[0][3][1] += R[0][3][2];
  446.     R[0][3][2] += R[0][3][3];
  447.  
  448.     R[1][0][0] += R[1][0][1];
  449.     R[1][0][1] += R[1][0][2];
  450.     R[1][0][2] += R[1][0][3];
  451.  
  452.     R[1][1][0] += R[1][1][1];
  453.     R[1][1][1] += R[1][1][2];
  454.     R[1][1][2] += R[1][1][3];
  455.  
  456.     R[1][2][0] += R[1][2][1];
  457.     R[1][2][1] += R[1][2][2];
  458.     R[1][2][2] += R[1][2][3];
  459.  
  460.     R[1][3][0] += R[1][3][1];
  461.     R[1][3][1] += R[1][3][2];
  462.     R[1][3][2] += R[1][3][3];
  463.  
  464.     R[2][0][0] += R[2][0][1];
  465.     R[2][0][1] += R[2][0][2];
  466.     R[2][0][2] += R[2][0][3];
  467.  
  468.     R[2][1][0] += R[2][1][1];
  469.     R[2][1][1] += R[2][1][2];
  470.     R[2][1][2] += R[2][1][3];
  471.  
  472.     R[2][2][0] += R[2][2][1];
  473.     R[2][2][1] += R[2][2][2];
  474.     R[2][2][2] += R[2][2][3];
  475.  
  476.     R[2][3][0] += R[2][3][1];
  477.     R[2][3][1] += R[2][3][2];
  478.     R[2][3][2] += R[2][3][3];
  479.  
  480.     R[3][0][0] += R[3][0][1];
  481.     R[3][0][1] += R[3][0][2];
  482.     R[3][0][2] += R[3][0][3];
  483.  
  484.     R[3][1][0] += R[3][1][1];
  485.     R[3][1][1] += R[3][1][2];
  486.     R[3][1][2] += R[3][1][3];
  487.  
  488.     R[3][2][0] += R[3][2][1];
  489.     R[3][2][1] += R[3][2][2];
  490.     R[3][2][2] += R[3][2][3];
  491.  
  492.     R[3][3][0] += R[3][3][1];
  493.     R[3][3][1] += R[3][3][2];
  494.     R[3][3][2] += R[3][3][3];
  495.     }
  496. }
  497.  
  498. /* ------------------------------------------------------------------------ */
  499.  
  500. /*
  501.  * Extract the k'th column of the tensor a into the matrix b.
  502.  */
  503. static void extract(
  504.   register Matrix b,
  505.   register Tensor a,
  506.   register int k)
  507. {
  508. b[0][0] = a[0][0][k];
  509. b[0][1] = a[1][0][k];
  510. b[0][2] = a[2][0][k];
  511. b[0][3] = a[3][0][k];
  512.  
  513. b[1][0] = a[0][1][k];
  514. b[1][1] = a[1][1][k];
  515. b[1][2] = a[2][1][k];
  516. b[1][3] = a[3][1][k];
  517.  
  518. b[2][0] = a[0][2][k];
  519. b[2][1] = a[1][2][k];
  520. b[2][2] = a[2][2][k];
  521. b[2][3] = a[3][2][k];
  522.  
  523. b[3][0] = a[0][3][k];
  524. b[3][1] = a[1][3][k];
  525. b[3][2] = a[2][3][k];
  526. b[3][3] = a[3][3][k];
  527. }
  528.  
  529. /* ------------------------------------------------------------------------ */
  530.  
  531. static float addemup(Matrix m)
  532. {
  533. return (m[0][0] + m[0][1] + m[0][2] + m[0][3] +
  534. m[1][0] + m[1][1] + m[1][2] + m[1][3] +
  535. m[2][0] + m[2][1] + m[2][2] + m[2][3] +
  536. m[3][0] + m[3][1] + m[3][2] + m[3][3]);
  537. }
  538.  
  539. /* ------------------------------------------------------------------------ */
  540.  
  541.